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Abstract 

In  this  study,  we  present  a  rigorous  3-D  mathematical  model,  to  treat  prediction  and  analysis  of  proton  exchange  membrane  fuel  cells 
(PEMFC)  species  concentration  and  current  density  distributions  in  different  flow  field  patterns  and  operating  conditions.  The  model  is  based 
on  the  solution  of  the  conservation  equations  of  mass,  momentum,  species  and  electric  current  in  a  fully  integrated  finite-volume  solver 
using  the  CFDRC  commercial  code.  The  polarization  curve  of  serpentine  flow  pattern  is  well  correlated  with  experimental  data.  The  cell 
performance  with  parallel  straight,  serpentine  and  interdigitated  flow  patterns  are  calculated  and  compared.  The  simulation  results  reveal  that 
serpentine  and  interdigitated  flow  patterns  show  strong  convection  and  high  mass  transfer.  However,  they  also  have  larger  pressure  loss.  In 
addition,  the  effects  of  operating  temperature  and  relative  humidity  are  also  studied.  Non-uniform  distributions  of  concentration  and  current 
density  appear  at  high  temperature,  high  current  density  and  low  humidity  operation,  which  could  lead  to  an  unstable  cell  performance. 

©  2005  Elsevier  B.V.  All  rights  reserved. 
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1.  Introduction 

Fuel  cell,  which  transforms  the  chemical  energy  into  elec¬ 
tricity  at  high  efficiency,  high  power  density,  is  considered 
to  be  a  promising  power  source,  especially  for  transporta¬ 
tion  application.  Its  performance  and  efficiency  must  be 
improved,  and  the  issues  of  cost,  reliability  and  safety  need  to 
be  considered  to  realize  the  fuel  cell  commerciality.  In  order 
to  enhance  its  performance  and  reliability,  it  is  necessary  to 
learn  more  about  the  mechanism  that  causes  the  performance 
loss,  such  as,  non-uniform  concentration,  current  density  dis¬ 
tributions,  high  ionic  resistance  due  to  dry  membrane,  or  high 
diffusive  resistance  due  to  the  flooding  on  the  cathode.  The 
flow  field  and  water/thermal  management  of  fuel  cell  need 
optimal  design  to  achieve  high  performance  and  reliability. 
The  numerical  modeling  and  dedicated  experimental  tech¬ 
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nique  development  will  be  the  effective  tools  to  improve  the 
optimal  design  of  fuel  cell  system. 

To  understand  the  mechanism  of  fluid  mechanics  and 
electrochemical  processes  within  a  fuel  cell,  the  first  two 
models  were  published  in  the  early  1990s  by  Springer  et 
al.  [1]  and  Bernardi  and  Verbrugge  [2].  These  models  are 
one-dimensional  and  only  account  for  simple  diffusive  mass 
transport  and  electrochemical  kinetics.  After  that,  several 
two-dimensional  models  have  been  presented  by  researchers 
[3-5].  Most  of  these  models  compute  the  flow  field  along 
a  single  channel  to  study  the  reaction  species  and  current 
density  distributions.  Results  of  polarization  curve  are  well 
correlated  with  experimental  data.  Um  et  al.  [6]  considered 
the  two-phase  flow  model  within  porous  diffusion  layer  and 
flow  channel  to  study  the  flooding  of  GDL  and  gas  channel, 
which  impedes  gas  flow  and  also  reduces  efficiency. 

Recently,  three-dimensional  models  based  on  die  com¬ 
putational  fluid  dynamics  approach  were  presented  and 
computed  by  commercial  software  [7,8]  or  self-developed 
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Nomenclature 

F  Faraday  constant 

z'f  in  porous  media,  the  current  flowing  through 

the  solid  parts  of  the  porous  matrix  (A) 
is  in  porous  media,  the  current  flowing  through 

the  pores  (A) 

jo  reference  current  (A  m-3) 

Ji  diffusion  flux 

U  velocity  (ms-1) 

Yi  mass-fractions  of  ith  species 

Greek  letters 

c^a  anode  kinetic  constant 

ac  cathode  kinetic  constant 

s  porosity 

k  permeability  (m2) 

[A]  the  near- wall  molar  concentration  of  the  react¬ 

ing  species 

li  viscosity  (m2  s-1) 

00 i  production  rates  of  ith  species  in  the  gas  phase 

p  pressure  (N  m-2) 

@s  solid  potential  (V) 

0F  fluid  potential  (V) 

r  shear  force  tension 


numerical  model  of  the  SIMPLE  algorithm  [9-11].  The 
three-dimensional  models  account  for  the  effect  of  the 
complex  geometry,  specifically  interdigitated  flow  field. 
They  allow  a  parametric  study  for  a  realistic  flow  field,  con¬ 
centration  and  current  distributions.  The  simulation  results 
are  well  compared  with  the  experimental  data  of  polarization 
curve.  However,  the  influence  of  flow  field  design  upon  con¬ 
centration  and  current  density  distributions  is  less  discussed. 

Mench  et  al.  [12]  proved  that  the  effects  of  cathode  stoi¬ 
chiometry  variation  and  transient  flooding  on  local  current 


density  affect  the  current  distribution.  When  voltages  are 
0.9-0.6V  (electrochemical  reaction  controlled  region),  the 
current  density  in  downstream  is  higher  than  upstream.  On 
the  other  hand,  when  the  voltages  are  0.4-0. 2  V,  the  current 
density  in  upstream  is  higher  than  downstream.  There  are 
significant  progresses  on  developed  experimental  technique 
for  the  measurement  of  current,  species,  temperature  distri¬ 
butions  [13]  and  flow  field  flooding  (water  content)  [14].  The 
efforts  on  fuel  cell  modeling  and  experimental  measurement 
technique  are  valuable  for  fuel  cell  developers,  which  can 
optimize  fuel  cell  designs  and  operations. 

The  major  objectives  of  this  work  are  to:  (1)  modify  a  3-D 
PEM  fuel  cell  model  to  predict  the  cell  performance,  con¬ 
centration  and  current  density  distributions  without  counting 
the  thermal  and  two-phase  flow  effect,  (2)  study  the  effect 
of  the  flow  field  patterns  on  concentration  and  current  den¬ 
sity  distributions,  which  can  aid  in  the  design  of  the  flow  field 
plates  and  (3)  study  the  cell  temperature  and  relative  humidity 
effects  on  the  uniformity  of  species  concentration  and  current 
density. 


2.  Model  development 

The  PEMFC  structure  includes  the  flow  field  plate,  gas 
diffusion,  catalyst  layer,  membranes  and  bipolar  plate,  etc. 
(see  Fig.  1).  Parallel  straight,  serpentine  flow  and  interdigi¬ 
tated  flow  patterns  are  shown  in  Fig.  2.  Tables  1  and  2  are 
the  dimensions  and  properties  of  the  flow  field  plate,  mem¬ 
brane,  electrode  material  properties  and  the  initial  operating 
conditions  used  in  the  numerical  simulation. 

The  numerical  calculations  are  based  on  solution  of 
the  conservation  equations  of  mass,  momentum,  energy, 
current  and  species  transport  on  a  computational  grid 
using  finite- volume,  finite-difference  methodologies.  The 
key  elements  to  modeling  fuel  cells  are  the  transport  phe¬ 
nomena  through  porous  media,  heterogeneous  reactions 


Current  collector 


Fig.  1.  Schematic  illustration  of  PEMFC. 


548 


F.-B.  Weng  et  al.  /  Journal  of  Power  Sources  145  (2005)  546-554 


Table  1 

Dimensions,  properties  and  parameters  for  the  numerical  model 


Channel  length  (mm) 

50 

Channel  width  (mm) 

1.2 

Channel  depth  (mm) 

1 

Rib  width  (mm) 

1.2 

Diffusion  layer  thickness  (mm) 

0.2 

Catalyst  layer  thickness  (mm) 

Anode:  0.018;  cathode:  0.026 

Membrane  thickness  (mm) 

0.035 

Total  reaction  area  (cm2) 

25 

Effective  diffusivity 

Bruggeman  model,  t  =  7  for 
membrane  Bruggeman  model, 
r  =  1.5  for  diffusion  and  cata¬ 
lyst  layer 

Membrane  permeability  (m2) 

1.8  x  10-18 

Diffusion  and  catalyst  layer  permeability 
(m2) 

1.76  x  10-11 

Membrane  porosity 

0.28 

Diffusion  and  catalyst  layer  porosity 

0.4 

Air  and  fuel  side  pressure  (atm) 

1 

Transfer  coefficient  (Tafel  constants)  at 
anode 

0.5 

Concentration  dependence  at  anode 

0.5(H2) 

Reference  current  density  at  anode 
((Am-3)(m3  (kg  mol  IT)-1)172) 

9.23  x  108 

Transfer  coefficients  (Tafel  constants)  at 
cathode 

1.5 

Concentration  dependence  at  cathode 

l-0(O2) 

Reference  current  density  at  cathode 
((Am-3)(m3  (kg  mol  IT)-1))172) 

1.05  x  106 

Diffusion  and  catalyst  layer  conductivity 
(l(Qm)-1) 

100 

within  porous  electrodes  and  the  coupling  between  mass 
transports,  electrochemical  reactions  and  current-potential 
fields. 


Table  2 

Operating  conditions  for  the  numerical  model 


H2  at  fuel  inlet  (cm3  min-1) 

300 

Anode  gas 

Humidified  hydrogen  (0,  50,  100%RH) 

02  at  fuel  inlet  (cm3  min- 1 ) 

300 

Cathode  gas 

Humidified  oxygen  (0,  50,  100%RH) 

Operating  pressure  (atm) 

1 

Operating  temperature  (°C) 

30,  50,  70 

The  main  assumptions  of  the  modeling  are: 

1.  steady  state; 

2.  laminar  flow; 

3.  isothermal; 

4.  all  gases  are  treated  as  ideal  gas; 

5.  gas  phase  flow  in  catalyst  layer,  gas  diffusion  layer  and 
channel  flow  field,  without  the  effects  of  vaporization 
and  condensation; 

6.  the  Stefan-Max  well  equations  for  multi- species  diffu¬ 
sion; 

7.  Butler- Volmer  equation  is  used  to  describe  electrochem¬ 
ical  reactions  within  the  catalyst  layers; 

8.  Nerst-Planck  equation  is  used  for  the  transport  of  protons 
through  the  membrane; 

9.  Ohm’s  law  is  applied  across  whole  region  of  the  polar¬ 
ization  curve; 

10.  gravity  effect  is  neglected. 

2.7.  Channel  flow  field 


According  to  assumptions,  these  conservation  equations 
of  mass,  momentum  and  species  in  flow  channel  are  follows: 


Mass  conservation: 
dp 


dt 


+  W(pU)  =  0 


Momentum  conservation: 
dpu 


dp 


A-axis  :  — - b  V  •  ( pUu )  =  — —  +  V(/z  Vw) 

dt  dx 

dpv  -»  dp 

Y- axis  : - j-  V  •  ( pUv )  = - b  V(/x  Vu) 

dt  dy 

dpw  dp 

Z- axis  : - b  V  •  ( pUw )  = - b  V(/z  Vw) 

dt  dz 


(1) 


(2) 


where  p  is  the  fluid  density,  p  the  pressure  and  /x  is  the 
dynamic  viscosity. 

Species  conservation: 


^(spYi)  +  V-(spUYi)  =  WJi 
ot 


(3) 


where  7/  is  the  mass-fraction  of  ith  species  and  //  is  the 
diffusive  flux. 
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The  species  diffusion  flux  may  be  written  as: 


Ji  —  pDi  VYi  + 


pYf 

M 


Dt  VM  -  pY/ffD,  VYj 

j 


AM 

~M~ 


where  M  is  the  mixture  molecular  weight  and  D[  is  the  effec¬ 
tive  mass  diffusion  coefficient  of  species  i. 


where  A,fs  is  the  free  stream  diffusion  coefficient  of  the  ith 
species,  s  porosity  of  the  medium  and  r  is  the  tortuosity  of 
the  medium.  It  is  a  common  practice  to  use  tortuosity  value  of 
1.5  in  Eq.  (5),  resulting  in  the  so-called  Bruggeman  model, 
and  the  same  value  was  adopted  in  our  simulations  due  to 
lack  of  better  information. 


2.2.  Gas  diffusion  layer 


The  main  assumption  is  isotropic  porous  media,  their  mass 
and  momentum  equations  are  effected  by  porosity  s  and  per¬ 
meability  K. 

These  equations  are  showed  as: 


Mass  conservation: 

■^-(sp)  +  V  •  (spU)  =  0 

ot 

where  p  is  the  fluid  density. 

Momentum  conservation: 

d  -  - 

—  ( [sp )  T  V  •  ( spU  •  t/)  —  — £V p  T  V  •  (£t)  T  sB  T 

dt 


s2fiU 

k 


where  p  is  the  fluid  density,  p  the  pressure,  fi  the  dynamic 
viscosity,  r  shear  force  tensor,  U  the  flow  velocity  vector,  s 
porosity  of  the  medium,  B  the  body  force  vector  and  k  is  a 
quantity  representing  the  square  of  the  effective  volume  to 
surface  area  ratio  of  the  porous  matrix. 

The  last  term  in  Eq.  (7)  represents  Darcy’s  drag  force 
imposed  by  pore  walls  on  the  fluid  within  the  pores,  and  usu¬ 
ally  results  in  a  significant  pressure  drop  across  the  porous 
medium.  Species  conservation  equation  in  the  gas  diffusion 
layer  is  the  same  as  in  flow  field.  The  continuity  of  current 
within  any  material  under  electroneutral  conditions  leads  to: 

V  •  i  =  0  (8) 

where  i  is  the  current  density  vector. 


equations  for  individual  gas  phase  species  (i.e.,  for  i  =  1,  . . ., 
Nq)  may  be  written  as: 

^(spYi)  +  V  •  (spUYi)  =  WJt  +  6>i  (9) 

ot 

where  F;  are  the  mass  fractions  of  the  ith  species  and  cot 
are  the  production  rates  of  the  ith  species  in  the  gas  phase. 
d>i  due  to  heterogeneous  electrochemical  reactions  within 
porous  media. 

The  species  diffusion  flux  may  be  written  as: 


pYi 


h  =  pD;  VK,  +  D,  VM  -  pY/ffDj  VYj 


M 

-PYi  — 


J 


M 


(10) 


j 


where  the  M  is  the  mixture  molecular  weight.  The  first  term 
represents  Fickian  diffusion  due  to  concentration  gradients. 
The  last  three  terms  are  correction  terms  necessary  to  satisfy 
the  Stefan-Maxwell  equations  for  multicomponent  species 
i  within  the  porous  medium,  and  depend  on  the  porosity,  s, 
and  tortuosity,  r,  of  the  medium: 


(11) 


where  D[  is  the  free  stream  diffusion  coefficient  of  the  ith 
species. 

In  the  case  of  electrochemical  reaction,  the  volumetric 
production  rate  of  a  given  species  is  expressed  by  the  ratio  of 
the  transfer  current,  jj,  and  the  Faraday  constant,  F: 


C0i  =  (a'!  -  a'iff  (12) 

r 

where  a"  and  a't  are  the  stoichiometric  coefficients  of  the 
products  and  reactants,  respectively.  The  transfer  current,  jj, 
is  obtained  from  the  Butler- Volmer  condition,  and  may  be 
written  in  its  most  general  form  as: 


ntl^.ref]0^ 


exp 


(  «a  JF 

V  RT 


—  exp 


RT 


N 

xj] \[AkY*-J 

k= 1 


(13) 


where  jo  is  the  reference  current  at  a  known  open  circuit  volt¬ 
age  and  a  &  j  and  otc  j  are  kinetic  constants.  [A]  represents  the 
near- wall  molar  concentration  of  the  reacting  species,  N  the 
total  number  of  reacting  species  and  a^j  are  the  concentration 
exponents.  Species  conservation  equation  in  the  gas  diffusion 
layer  is  the  same  in  flow  field.  The  continuity  of  current  within 
any  material  under  electroneutral  conditions  leads  to: 


V  •  i  =  0 


(14) 


2.3.  Catalyst  layer 

The  mass  and  momentum  equation  in  the  catalyst  layer  is 
the  same  in  the  gas  diffusion  layer.  The  mass  conservation 


where  i  is  the  current  density  vector.  When  the  material  is 
a  porous  electrode,  the  current  may  split  into  two  parts:  one 
flowing  through  the  polymer  electrolyte  (ionic  phase)  and  the 
other  flowing  through  the  solid  parts  or  electronic  phase  of 
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the  porous  matrix.  If  these  components  are  denoted  by  /p  and 
is,  respectively,  Eq.  (14)  can  be  rewritten  as: 

V/F  +  Vis  =  0  (15) 


During  electrochemical  reactions  within  a  porous  solid, 
flow  of  current  results  from  charge  transfer  between  the  ionic 
and  electronic  phase.  This  results  in  a  source  in  one  of  the 
phases  and  a  sink  in  the  other  phase.  Thus,  Eq.  (15)  may  be 
written  as: 


-V  •  z'f  +  V  •  is  =  j't 


(16) 


[S/  V]eff  is  the  effective  surface-to- volume  ratio. 
Application  of  Ohm’s  law  to  Eq.  (15),  yields: 


V  •  (crp  V 0p)  = 


-V  •  (as  V^s)  =  Jt 


(17) 


where  <£>f  and  0s  are  the  electric  potentials  of  the  ionic  and 
electronic  phases,  respectively. 


2.4.  Proton  exchange  membrane 

In  the  proton  exchange  membrane,  governing  equations 
are  the  same  as  in  the  catalyst  layer,  but  without  any  chemical 
reaction.  The  membrane  proton  conductivity  is  correlated  by 
the  experimental  data,  Springer  et  al.  [  1] ,  which  are  dependent 
on  temperature  and  water  content. 

The  solution  to  the  above  set  of  conservation  equations 
is  performed  using  a  finite-volume  scheme  on  arbitrary 
unstructured  mesh  topology  within  the  framework  of  the 
commercial  CFD  code  CFD-ACE+,  provided  by  CFDRC 
software. 


2.5.  Computational  domain 

We  describe  a  detailed  three-dimensional  model  of  trans¬ 
port  phenomena  and  performance  within  the  PEMFC  includ¬ 
ing  gas  channel,  gas  diffusion  layer,  catalyst  layer  and 
membrane.  The  dimensions  of  the  computational  domains 
have  27  elements  in  the  v-direction,  84  elements  in  the 
y-direction,  42  elements  in  the  z-direction,  for  a  total  of 
about  95,256  elements  and  198,660  nodes.  The  gird  dis¬ 
tributions  and  cross-section  of  concentration  and  current 
density  are  shown  in  Su  et  al.  [8].  The  species  concentra¬ 
tions  are  calculated  at  the  interface  of  catalyst/GDL,  while 
the  current  density  distributions  are  at  the  current  collect 
plate. 

2.6.  Experimental  details 

An  experimental  system  is  setup  to  obtain  the  cur¬ 
rent/voltage  polarization  curves.  Graphite  plate  is  used  as 
a  current  collector  with  serpentine  flow  pattern.  Other  mate¬ 
rials  used  in  the  single  cell  are  Gore  5621  ME  A,  reaction 
area  25  cm2  (from  Gore  Tex  Fuel  Cell  Co.)  and  brass  current 


Current  Density  (mA  /cm^) 


Fig.  3.  Comparison  of  experiment  with  simulation  polarization  curves  (ser¬ 
pentine  flow  pattern). 

collector.  The  test  cell  configuration  is  shown  in  Fig.  1.  The 
dimensions  of  components  and  operating  conditions  are  con¬ 
sistent  with  the  numerical  model.  The  fuel  cell  test  system  is 
developed  by  Electrochem.  Co. 


3.  Results  and  discussions 

The  results  of  numerical  simulation  are  validated  with 
experimental  data.  The  polarization  curves  of  serpentine  flow 
pattern  at  base  operating  condition  are  calculated,  as  shown 
in  Fig.  3.  The  properties,  parameters  and  base  conditions  of 
numerical  simulation  are  listed  in  Tables  1  and  2,  which  are 
referred  mainly  from  Bernardi  and  Verbrugge  [2].  The  rel¬ 
ative  deviation  of  current  density  of  the  modeling-predicted 
data  with  experimental  data  is  less  than  20%  under  high  volt¬ 
age  0.7  V,  and  5%  for  low  voltage  0.4  V.  Since  we  neglect 
water  flooding  effect  at  the  cathode  side  in  the  single-phase 
model,  some  data  are  slightly  over  predicted.  However,  the 
trend  of  increasing  current  density  with  decreasing  the  cell 
potential  agrees  well. 

The  polarization  curves  with  various  flow  field  patterns 
are  shown  in  Fig.  4.  The  design  of  parallel  straight  pat¬ 
tern  has  less  power  density,  around  2/3,  comparing  with 
serpentine  and  interdigitated  flow  patterns.  In  addition,  the 
interdigitated  flow  pattern  can  enhance  limiting  current  den¬ 
sity  due  to  the  convective  flow  field  at  porous  diffusion 
layer. 


Fig.  4.  Polarization  curves  of  straight  flow,  serpentine  flow  and  interdigitated 
flow  patterns  in  the  same  operating  conditions  (numerical  simulation). 
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Fig.  5.  Parallel  straight  flow  (a)  O2  mass  fraction,  (b)  current  density,  serpentine  flow,  (c)  H2  mass  fraction,  (d)  current  density  and  interdigitated  flow,  (e)  H2 
mass  fraction  and  (f)  current  density  (0.6  V). 


Fig.  5  displays  the  species  concentrations  and  the  current 
density  distributions  (operated  at  0.6  V).  Fig.  5(a)  illustrates 
the  non-uniform  oxygen  concentration  distribution  at  cathode 
side  for  the  parallel  straight  pattern.  Lower  oxygen  concen¬ 


tration  appears  at  cell  center  and  left-top  corner  resulting  low 
current  density.  Fig.  5(c  and  d)  displays  hydrogen  concen¬ 
tration  and  current  density  of  the  serpentine  flow  pattern. 
The  downstream  of  channel  has  lower  hydrogen  concentra- 
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tion,  which  also  resulting  low  current  density.  Fig.  5(e  and 
f)  displays  the  hydrogen  concentration  and  current  density 
for  interdigitated  flow  pattern.  The  hydrogen  concentration 
and  current  density  are  higher  and  more  uniform  at  inlet 
channel.  On  the  contrary,  the  concentrations  are  significantly 
decreased  at  porous  diffusion  layer.  The  hydrogen  concen¬ 
tration  is  depleted  near  outlet  corner.  The  major  pressure  loss 
also  appears  at  the  cross-section  of  porous  diffusion  layer. 
It  is  noticed  that  the  decrease  of  the  rib  width  of  flow  chan¬ 
nel  will  cause  the  decrease  of  the  pressure  and  concentration 
drop.  This  leads  to  the  uniformity  of  the  concentration  as 
well  as  the  increase  of  cell  performance.  We  suggest  that  the 
optimal  value  of  the  rib  width  over  diffusion  layer  thickness 
around  two.  As  the  operating  potential  is  smaller  than  0.4  V, 
the  non-uniform  effect  becomes  more  pronounced.  Even¬ 
tually,  the  membrane  could  dry  out  at  high  current  density 
region. 

Fig.  6  shows  the  numerical  prediction  of  polarization 
curves  at  variable  operating  cell  temperatures.  The  cell  per¬ 
formance  increases  with  temperature  due  to  high  membrane 
conductivity.  This  is  contradicted  with  experimental  results 


Fig.  6.  Polarization  curves  of  the  fuel  cell  at  variable  operating  temperatures 
(numerical  simulation). 

[10]  with  operating  temperature  above  70  °C.  In  the  model, 
the  membrane  dry  out  situation  is  not  considered  at  high 
cell  temperature.  However,  higher  fuel  and  oxygen  tem¬ 
peratures  could  also  induce  more  water  vapor,  resulting  in 
lower  hydrogen  concentration  at  inlet.  Therefore,  the  lim¬ 
iting  current  density  decreases  as  fuel/oxygen  temperatures 
increase. 


Fig.  7.  H2  mass  fraction  at  operation  voltage  (0.2  V)  and  temperature:  (a)  30  °C  and  (b)  70  °C. 
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Fig.  8.  Current  density  at  operation  voltage  0.2  V  and  temperature:  (a)  30  °C  and  (b)  70  °C. 
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Figs.  7  and  8  display  the  hydrogen  concentration  and 
current  density  at  30  and  70  °C,  respectively.  The  operating 
conditions  are  0.2  V,  1  atm  and  flow  rate  is  300  cm3  min-1. 
The  high  temperature  of  fuel  cell  has  low  ohmic  resistance, 
and  more  active  electrochemical  reaction  at  upstream  of 
flow  channel.  Therefore,  the  limiting  current  appears  at 
downstream  of  flow  channel.  However,  the  concentration 
and  current  density  become  non-uniform.  The  mem¬ 
brane  could  dry  out  at  high  current  density  area  due  to 
the  heat  generation,  which  could  lead  an  unstable  cell 
performance. 

Fig.  9  shows  the  polarization  curves  of  the  fuel  cell  at 
variable  relative  humidity.  As  the  relative  humidity  is  above 
50%,  the  membrane  conductivity  is  constant.  The  100%  rel¬ 
ative  humidity  has  higher  limiting  current  than  that  of  50% 
humidity.  The  fuel  cell  needs  high  relative  humidity  as  the  fuel 
cell  operating  high  current  density.  The  membrane  conduc¬ 
tivity  is  a  function  of  water  content  and  temperature,  which 
proposed  by  Springer  et  al.  [1].  The  numerical  results  show 
that  membrane  is  dry  out  and  high  ohmic  resistance  without 
humidified  fuel/oxygen  gases. 

Figs.  10  and  11  display  the  hydrogen  concentration  and 
current  density  distributions  at  low  voltage  (0.2  V)  and  vari- 
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Fig.  9.  Polarization  curves  of  the  fuel  cell  at  variable  relative  humidity 
(numerical  simulation). 


able  relative  humidity.  Without  humidified  hydrogen/oxygen 
gases,  the  membrane  has  low  conductivity  and  low  cur¬ 
rent  density.  The  current  density  at  the  downstream  of  flow 
channel  increases  due  to  the  water  generation  by  chemical 
reaction.  On  the  other  hand,  with  100%  relative  humidity  of 
hydrogen/oxygen  gases,  the  hydrogen  is  fast  consumed  and 
depleted  at  downstream  of  flow  channel.  This  leads  to  the 
decrease  of  current  density.  The  operating  condition  of  low 


Fig.  10.  H2  mass  fraction  at  operation  voltage  0.2  V:  (a)  RH  =  0%  and  (b)  RH  =  100%. 
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Fig.  11.  Current  density  at  operation  voltage  0.2  V  (a)  RH  =  0%  and  (b)  RH  =  100%. 
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humidity  and  high  current  density  could  lead  an  unstable  cell 
performance. 

4.  Conclusions 

A  3-D  mathematical  model  to  predict  species  concentra¬ 
tion,  current  density  distributions  as  well  as  cell  performance 
is  presented  in  this  work.  The  model  is  based  on  the  solution 
of  the  conservation  equations  of  mass,  momentum,  species 
and  electric  current  in  a  fully  integrated  finite-volume  solver 
using  the  CFDRC  commercial  code.  The  influence  of  flow 
field  pattern,  cell  temperature  and  humidity  on  the  non¬ 
uniformity  of  species  concentration,  current  is  investigated. 
The  polarization  curve  of  serpentine  flow  pattern  is  well  cor¬ 
related  with  experimental  data.  The  following  conclusions  are 
made: 

1 .  The  simulation  results  reveal  that  serpentine  and  interdigi- 
tated  flow  patterns  show  strong  convection  and  high  mass 
transfer.  However,  they  also  have  higher  pressure  loss. 
In  the  interdigitated  flow  pattern,  the  cell  performance, 
non-uniformity  and  pressure  drop  are  strong  dependent 
on  the  width  of  channel  rib.  Optimal  design  is  important 
to  improve  its  cell  performance. 

2.  The  cell  performance  increases  with  the  temperature 
increase  due  to  the  high  membrane  conductivity.  How¬ 
ever,  the  membrane  will  dry  out  at  high  cell  temperature. 

3.  Non-uniform  distributions  of  gas  concentration  and  cur¬ 
rent  density  appear  at  high  temperature,  high  current  den¬ 
sity  and  low  humidity  operations,  which  could  lead  to  an 
unstable  cell  performance. 
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